Projekt zakładał analizę danych pacjentów chorych na COVID-19. Przetworzenie danych ze zbioru umożliwiło sporządzenie raportu zawierającego wykresy takie jak rozkład wyzdrowień i śmierci. Podjęto także próbę stworzenia klasyfikatora, który będzie potrafił wskazać czy dany pacjent przeżyje.
library("readxl")
library("dplyr")
library(skimr)
library(plotly)
library(ggplot2)
library(DT)
library(hablar)
library(reshape2)
library(tidyr)
library(caret)
set.seed(23)
Data próbki została przycięte do formatu “yyyy-MM-dd”.Następnie dane zostały zgrupowane pod dacie przyjęcia oraz dacie opuszenia szpitala, z powodu wyzdrowienia lub zgonu. Co umożliwiło przypisanie brakujących PATIENT_ID.
df <- mutate(data, RE_DATE = as.Date(substr(RE_DATE,1,10), "%Y-%m-%d"))
df <- df %>% group_by(Admission.time, Discharge.time) %>% mutate(PATIENT_ID = max(PATIENT_ID, na.rm = TRUE)) %>% ungroup()
df <- df %>% group_by(PATIENT_ID, RE_DATE) %>%
mutate(Admission.time = max(Admission.time, na.rm = TRUE)) %>%
mutate(Discharge.time = max(Discharge.time, na.rm = TRUE))%>% ungroup()
Wszystkie próbki z danego dnia dla danego pacjenta zostały zgrupowane.
df <- df %>% group_by(PATIENT_ID, RE_DATE, Admission.time, Discharge.time, age, outcome, gender) %>% summarise_all(sum_) %>% ungroup()
Dane liczbowe odpowiadające płci zostały zamienione na dane tekstowe.
df <- df %>% mutate(gender = factor(gender, labels = c('male', 'female')))
df <- df %>% mutate(outcome = factor(outcome, labels = c('survived', 'died')))
Kolumny i wiersze zawierające niekompletne dane zostały usunięte.
df$X2019.nCoV.nucleic.acid.detection = NULL
df <- subset(df, !is.na(RE_DATE))
Tabela przedstawia zbiór po przetworzeniu.
datatable(df,options = list(scrollX = TRUE), filter = 'top')
my_skim <- skim_with(character = sfl(), append = FALSE, base = sfl(n_missing))
datatable(my_skim(df), options = list(scrollX = TRUE))
Wykresy przedstawiają wpływ wartości dehydrogenazy mleczanowej (LDH), białka C-reaktywnego (hs-CRP), oraz liczby limfocytów.
last_data <- df
last_data <- last_data %>% group_by(PATIENT_ID) %>%
summarise(outcome = last(outcome), LDH = last(Lactate.dehydrogenase), Limfocyty = last(Lactate.dehydrogenase), 'hs-CRP' = last(High.sensitivity.C.reactive.protein)) %>% ungroup()
last_data[ ,c('PATIENT_ID')] <- list(NULL)
datatable(last_data)
last_plot_data = melt(last_data, id=c("outcome"))
graph <- ggplot(last_plot_data, aes(x=value, y=outcome)) + geom_point() + facet_wrap(~variable)
ggplotly(graph)
Wykres przedstawia dane dla każdego dnia. Zsumowana została liczba pacjentów opuszczających szpital z podziałem na zgony oraz wyzdrowienia.
by_date <- mutate(df, Discharge.time = substr(Discharge.time,1,10))
by_date <- by_date %>% group_by(Discharge.time, PATIENT_ID) %>% mutate(outcome = first(outcome == "died")) %>%
summarise(outcome = first(outcome)) %>% group_by(Discharge.time) %>%
summarise(Pacjenci = n_distinct(PATIENT_ID), Zgony = sum(outcome == TRUE), Wyzdrowienia = sum(outcome == FALSE))
plot_data = melt(by_date, id=c("Discharge.time"))
p <- ggplot(plot_data) +
geom_point(aes(x = Discharge.time, y = value, color = variable,
text = paste(' Dzień: ', Discharge.time,'<br>','Liczba: ',value,'<br>','Typ: ',variable))) +
scale_colour_manual(values=c("black","red","green")) +
labs(color = "Legenda", x = "Dzień", y = "Liczba") +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1))
g <- ggplotly(p, tooltip = "text") %>% layout(legend = list(x = 0.05, y = 0.95))
g
Przed przystąpienie do uczenia klasyfikatora dane zostały odpowiednio przygotowane. Ze zbioru zostali usunięci pacjenci ze zbyt mała liczbą przeprowadzonych badan oraz badańia wykonane na zbyt małej liczbie pacjentów.
raw_train_data <- df
raw_train_data[ ,c('Admission.time', 'Discharge.time')] <- list(NULL)
raw_train_data <- raw_train_data[which(rowMeans(!is.na(raw_train_data)) > 0.7), which(colMeans(!is.na(raw_train_data)) > 0.32)]
Po wstępnym przetworzeniu zbioru, usunięto wszystkich pacjentów, którzy nie mieli kompletu pozostałych badań.
raw_train_data <- raw_train_data %>% drop_na()
patient_data <- raw_train_data %>% group_by(PATIENT_ID) %>% summarise_at(vars(outcome, gender), last) %>% ungroup()
raw_train_data <- raw_train_data %>% group_by(PATIENT_ID) %>% summarise_if(is.numeric, max) %>% ungroup()
Tabela przedstawia ostateczny zbiór, który zostanie wykorzystany do uczenia klasyfikatora. Zawiera on 217 rekordóW.
raw_train_data <- merge(x = patient_data, y = raw_train_data, by ="PATIENT_ID", all = TRUE)
datatable(raw_train_data,options = list(scrollX = TRUE))
Dane treningowe stanowią 75% zbioru. Pozostałe 25% zostało podzielone pomiędzy dane testowe (40%) i walidujące (60%).
inTraining <- createDataPartition(y = raw_train_data$outcome, p = .75, list = FALSE)
training <- raw_train_data[inTraining,]
testing_val <- raw_train_data[-inTraining,]
testing_val_index <- createDataPartition(y = testing_val$outcome, p = .6, list = FALSE)
testing <- testing_val[-testing_val_index,]
val <- testing_val[testing_val_index,]
datatable(training, options = list(scrollX = TRUE))
datatable(testing, options = list(scrollX = TRUE))
datatable(val, options = list(scrollX = TRUE))
Dane zostały wykorzystane do uczenia klasyfikator random forest.
ctrl <- trainControl(method = "repeatedcv", number = 2,repeats = 5)
fit <- train(outcome ~ .,
data = training,
method = "rf",
trControl = ctrl,
ntree = 10)
fit
## Random Forest
##
## 163 samples
## 58 predictor
## 2 classes: 'survived', 'died'
##
## No pre-processing
## Resampling: Cross-Validated (2 fold, repeated 5 times)
## Summary of sample sizes: 81, 82, 81, 82, 82, 81, ...
## Resampling results across tuning parameters:
##
## mtry Accuracy Kappa
## 2 0.9164408 0.8305415
## 30 0.9839958 0.9676142
## 58 0.9950918 0.9901007
##
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was mtry = 58.
rfGrid <- expand.grid(mtry = 10:30)
gridCtrl <- trainControl(method = "repeatedcv", summaryFunction = twoClassSummary, classProbs = TRUE, number = 2, repeats = 5)
fitTune <- train(outcome ~ .,
data = training,
method = "rf",
metric = "ROC",
preProc = c("center", "scale"),
trControl = gridCtrl,
tuneGrid = rfGrid,
ntree = 30)
fitTune
## Random Forest
##
## 163 samples
## 58 predictor
## 2 classes: 'survived', 'died'
##
## Pre-processing: centered (58), scaled (58)
## Resampling: Cross-Validated (2 fold, repeated 5 times)
## Summary of sample sizes: 81, 82, 82, 81, 81, 82, ...
## Resampling results across tuning parameters:
##
## mtry ROC Sens Spec
## 10 0.9933192 0.9666667 0.9669670
## 11 0.9986728 0.9777778 0.9835586
## 12 0.9976852 0.9688889 0.9807057
## 13 0.9978095 0.9755556 0.9863363
## 14 0.9971330 0.9688889 0.9725976
## 15 0.9971296 0.9755556 0.9694444
## 16 0.9979630 0.9844444 0.9780030
## 17 0.9983025 0.9755556 0.9888889
## 18 0.9981790 0.9666667 1.0000000
## 19 0.9982716 0.9777778 0.9917417
## 20 0.9991358 0.9822222 0.9945195
## 21 0.9984259 0.9733333 0.9833333
## 22 0.9989506 0.9888889 0.9889640
## 23 0.9983959 0.9844444 0.9972973
## 24 0.9988889 0.9822222 0.9889640
## 25 0.9992009 0.9866667 1.0000000
## 26 0.9990741 0.9844444 0.9972222
## 27 0.9999074 0.9866667 1.0000000
## 28 0.9991975 0.9844444 1.0000000
## 29 0.9996914 0.9888889 0.9972973
## 30 0.9992593 0.9888889 1.0000000
##
## ROC was used to select the optimal model using the largest value.
## The final value used for the model was mtry = 27.
Uzyskano bardzo wysoki wskaźnik “accuracy”. Nauczony klasyfikator poprawnie sklasyfikował wszystkie próbki. Niestety liczba danych wykorzystana w procesie jest za mała. Nie można mieć pewności, że klasyfikator będzie działać z tak dobrą skutecznością na znacznie większym zbiorze.
predict_result <- predict(fitTune, newdata = val)
confusionMatrix(data = predict_result, val$outcome)
## Confusion Matrix and Statistics
##
## Reference
## Prediction survived died
## survived 18 0
## died 0 15
##
## Accuracy : 1
## 95% CI : (0.8942, 1)
## No Information Rate : 0.5455
## P-Value [Acc > NIR] : 2.056e-09
##
## Kappa : 1
##
## Mcnemar's Test P-Value : NA
##
## Sensitivity : 1.0000
## Specificity : 1.0000
## Pos Pred Value : 1.0000
## Neg Pred Value : 1.0000
## Prevalence : 0.5455
## Detection Rate : 0.5455
## Detection Prevalence : 0.5455
## Balanced Accuracy : 1.0000
##
## 'Positive' Class : survived
##